February 19, 2015

What can we do in R?

  1. Read/Write Images
  2. Visualization of Images
  3. Inhomogeneity/Bias Field Correction
  4. Skull Stripping/Brain Extraction
  5. Image Registration
  6. Tissue-Class Segmentation
  7. Image operations

Introduction: Packages

Medical Imaging Task View

  • oro.nifti: read/write data, the nifti object
  • fslr: process data (need FSL for most of the functionality)
  • ANTsR: process data (full toolbox)
  • extrantsr: makes ANTsR work with nifti objects
  • dti - adaptive smoothing and diffusion tensor tools

  • fmri - post-processing analysis: linear models and p-value smoothing
  • AnalyzeFMRI - fMRI analysis (last updated in 2013)

Data used from NITRC

Multi-modal dataset from HAMMER, (NIfTI conversion from ANALYZE).

files
            t1             t2             pd          flair 
   "T1.nii.gz"    "T2.nii.gz"    "PD.nii.gz" "FLAIR.nii.gz" 

Basics: Read in the Files!

fslr: readnii uses oro.nifti::readNIfTI:

library(fslr)
base_t1 = readnii(files["t1"])
  • like an array
  • ANTsR uses pointers (faster), but not as intuitive

Quick Plots

library(fslr)
fslr::ortho2(base_t1)

Quick Plots: Overlays

over_50 = mask_img(base_t1, base_t1 > 40)
fslr::ortho2(base_t1, over_50)

Visualization: Image Slices

image(base_t1, z = 55, plot.type = "single")

Quick Plots: Overlays

over_50[over_50 <= 0] = NA; over_50 = cal_img(over_50)
overlay(base_t1, over_50, z = 55, plot.type = "single")

Bias Field Correction

ANTsR/extrantsr

  • bias_correct from extrantsr package calls ANTsR::n4BiasFieldCorrection
  • EM-like, assumes bias is smooth over space, logs the data
library(extrantsr)
n4_t1 = bias_correct(file = base_t1, correction = "N4", retimg = TRUE)

fslr: Uses method by Sled, Zijdenbos, and Evans (1998) (slow)

bc_t1 = fsl_biascorrect(file = base_t1)
FSLDIR='/usr/local/fsl/'; export FSLDIR; sh "${FSLDIR}/etc/fslconf/fsl.sh"; FSLOUTPUTTYPE=NIFTI_GZ; export FSLOUTPUTTYPE; $FSLDIR/bin/fast    -B --nopve --out="/var/folders/1s/wrtqcpxn685_zk570bnx9_rr0000gr/T//Rtmpa4mKOM/file15e6f7f9ff9a7" "/private/var/folders/1s/wrtqcpxn685_zk570bnx9_rr0000gr/T/Rtmpa4mKOM/file15e6f60841581.nii.gz";

Bias Field Correction: Results

ratio = finite_img(n4_t1 / base_t1)
ortho2(n4_t1, ratio, col.y = alpha(hotmetal(), 0.5))

Skull Stripping: FSL BET

ss_t1 = fslbet(n4_t1, outfile = "SS_Image")
FSLDIR='/usr/local/fsl/'; export FSLDIR; sh "${FSLDIR}/etc/fslconf/fsl.sh"; FSLOUTPUTTYPE=NIFTI_GZ; export FSLOUTPUTTYPE; $FSLDIR/bin/bet2 "/private/var/folders/1s/wrtqcpxn685_zk570bnx9_rr0000gr/T/Rtmpa4mKOM/file15e6f2126f0db.nii.gz" "./SS_Image" 

Overlaying Skull Stripped mask

mask = ss_t1 > 0
ortho2(base_t1, y = mask, col.y = alpha("red", 0.5))

Visualization: Cropped Image Slices

cropped = dropEmptyImageDimensions(ss_t1)
image(cropped, z = floor(dim(cropped)[3]/2), plot.type = "single")

Visualization: 3-dimensions

devtools::source_gist("bd40d10afabc503d71e8")

Rigid-Body Registration

ANTsR/extrantsr

  • antsRegistration - rigid/affine/non-linear diffeomorphic
  • extrantsr::registration - wraps antsRegistration to use nifti objects

fslr

  • flirt - linear/affine registration
  • fnirt - non-linear registration (need affine first)
  • fnirt_with_affine - wraps above 2

Image Registration: Rigid-Body Registration

ants_reg_flair = registration(
  filename = files["flair"], 
  template.file = n4_t1, 
  typeofTransform = "Rigid")

Rigid-Body Registration Results

double_ortho(n4_t1, ants_reg_flair$outfile)

Rigid-Body Registration Slice

multi_overlay(list(n4_t1, ants_reg_flair$outfile), z = 55)

Non-linear Registration

template.file = mni_fname(mm = "1", brain = TRUE)
ss_t1_to_mni = registration(
  filename = ss_t1, 
  template.file = template.file, 
  typeofTransform = "SyN", remove.warp = FALSE,
  outprefix = "temp")

Non-linear Registration Results

double_ortho(template.file, ss_t1_to_mni$outfile)

Applying Registration Transformations

reg_flair_to_mni = ants_apply_transforms(
  fixed = template.file, 
  moving = ants_reg_flair$outfile, # registered FLAIR
  interpolator = "Linear",
  transformlist = ss_t1_to_mni$fwdtransforms )

Non-linear Registration Results FLAIR

double_ortho(template.file, reg_flair_to_mni)

Preprocess MRI: within a visit

extrantsr::preprocess_mri_within will do inhomogeneity correction, skull strip (or mask), and register to the first scan.

proc_images = preprocess_mri_within(
  files = files[c("t1", "t2", "pd", "flair")], 
  maskfile = ss_t1 > 0)

Preprocess MRI: across visits

preprocess_mri_across combines preprocess_mri_within and registration. If you had baseline/follow-up data:

outfiles = gsub("[.]nii", '_process.nii', files)
preprocess_mri_across(
  baseline_files = files[c("base_t1", "base_t2", "base_pd", "base_flair")],
  followup_files = files[c("f_t1", "f_t2", "f_pd", "f_flair")],
  baseline_outfiles = outfiles[c("base_t1", "base_t2", "base_pd", "base_flair")],
  followup_outfiles = outfiles[c("f_t1", "f_t2", "f_pd", "f_flair")],
  maskfile = "Brain_Mask.nii.gz")

Tissue-Class Segmentation

fslr: uses FAST

  • --nobias as an option does not do bias field correction (if already done)
fast_t1 = fast(ss_t1, opts = "--nobias")
FSLDIR='/usr/local/fsl/'; export FSLDIR; sh "${FSLDIR}/etc/fslconf/fsl.sh"; FSLOUTPUTTYPE=NIFTI_GZ; export FSLOUTPUTTYPE; $FSLDIR/bin/fast   --nobias --out="/var/folders/1s/wrtqcpxn685_zk570bnx9_rr0000gr/T//Rtmpa4mKOM/file15e6f60cbfee3" "/private/var/folders/1s/wrtqcpxn685_zk570bnx9_rr0000gr/T/Rtmpa4mKOM/file15e6f682913b1.nii.gz";

ANTsR/extrantsr: uses Atropos

  • antsRegistration - atropos
  • extrantsr - otropos
tissue_seg = otropos(
  a = ss_t1,
  x = mask)

Tissue-Class Segmentation Results

double_ortho(ss_t1, tissue_seg$segmentation)

Image operations

fslr

  • fslsmooth - Gaussian/box smoothing
  • fslerode/fsldilate - erosion/dilation
  • fslfill/fslfill2 - fill holes

spm12r

  • spm_bwlabel - label connected components

ANTsR

  • smooth_image - Gaussian smoothing
  • oMath("ME")/oMath("MD") - erosion/dilation
  • oMath("FillHoles") - fill holes
  • oMath("GetLargestComponent") - find largest components

fMRI

  • fsl_slicetimer - slice timing correction
  • ANTsR::preprocessfMRI
  • spm12r package calls out MATLAB using SPM

Intensity Normalization

  • WhiteStripe
  • Whole brain z-scoring
  • Histogram matching
  • General standardization methods

Overview

  • Methods are being developed for neuroimaging analysis in R
  • These are not standardized nor centralized
  • We need something like BioConductor
    • standard data structures
    • standard data sets
    • NITRC isn't exactly what we want

Combining Registrations

reg_flair_to_mni_1interp = ants_apply_transforms(
  fixed = template.file, 
  moving = files["flair"], # Raw FLAIR
  transformlist = c(ants_reg_flair$fwdtransforms,
                    ss_t1_to_mni$fwdtransforms),
  interpolator = "Linear")

Non-linear Registration

tfile = "FLAIR2.nii.gz"
reg_flair_to_mni_1interp_1step = registration(
  filename = ss_t1, 
  template.file = template.file, 
  typeofTransform = "SyN", 
  other.files = files['flair'],
  other.outfiles = tfile,
  other.init = ants_reg_flair$fwdtransforms)

Difference in Interpolations

rat2 = abs((reg_flair_to_mni_1interp - reg_flair_to_mni) / (reg_flair_to_mni_1interp +  reg_flair_to_mni)) 
rat2 = finite_img(rat2)
ortho2(template.file, rat2, col.y = alpha(hotmetal(), 0.5))

Difference in Interpolations Registrations

double_ortho(reg_flair_to_mni_1interp, reg_flair_to_mni)

Sled, John G, Alex P Zijdenbos, and Alan C Evans. 1998. “A Nonparametric Method for Automatic Correction of Intensity Nonuniformity in MRI Data.” Medical Imaging, IEEE Transactions on 17 (1). IEEE: 87–97.